Theory of nonlinear rheology and yielding of dense colloidal suspensions 
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A first principles approach to the nonlinear flow of dense suspensions is presented which captures 
shear thinning of colloidal fluids and dynamical yielding of colloidal glasses. The advection of density 
fluctuations plays a central role, suppressing the caging of particles and speeding up structural 
relaxation. A mode coupling approach is developed to explore these effects. 
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The properties of dispersions under flow are central to 
their processing and technological use |J 9. But espe- 
cially the non-linear rheology is not yet well understood. 
For the simplest case of steady shearing, the low density 
behavior is known || , but upon increasing the density the 
growing importance of particle interactions requires the- 
oretical approximation || , hinders simulations || , and 
calls for studies of model systems, e.g. fj], [§• Of major 
interest is the arrest of the structural relaxation when ap- 
proaching solidification for higher densities, which raises 
the question of how the imposition of steady shearing 
might interfere with glass formation. The linear phe- 
nomenology is familiar: a colloidal fluid possesses a vis- 
cosity and flows, while a colloidal glass characterized by 
elastic constants, only distorts under strain (l], |^j. But 
the nonlinear rheology of glassy colloids, which exhibit 
a continuous slowing down of the structural relaxation 
due to particle blocking (the "cage effect") jjj, is less 
clear. While the mode coupling theory (MCT) recovers 
the linear phenomenology of this fluid-to-glass transition 
from microscopic starting points a nonlinear exter- 
nal driving introduces new time scales whose influence 
on (non-)equilibration is not understood, and has been 
addressed only in minimal models [ flO| or mean-field ap- 
proaches Moreover, as the true nature of the glass is 
still uncertain, its behavior under shearing may provide 
broader new insights (as suggested by recent simulation 
studies £H| [if]). 

Here we develop a first-principles approach for the 
simplest case of a disordered colloidal suspension under 
steady imposed shear, neglecting both many-body hy- 
drodynamics and the resulting velocity fluctuations. We 
first identify some generic features in the yield properties 
of glass; approximations suggested by the MCT are then 
introduced in order to derive quantitative predictions. 

The system consists of N spherical particles (diameter 
d) dispersed in a volume V of solvent with imposed flow 
profile v(r) — nr, where for simple shear with velocity 
along the x-axis and its gradient along the y-axis, the 
shear rate tensor is Kij — 7 Si X 6j y . The effect of the 
shear rate 7 on the particle dynamics is measured by 
the Peclet number Pco = -jcP/Dq, formed with the 
bare diffusion coefficient Dq of a particle. Dimensionless 



units are obtained by setting d = Dq = ksT = 1. The 
evolution of the distribution function of the particle 
positions, r^, i = 1,...N, under internal forces F,; and 
shearing, but neglecting hydrodynamic interactions, is 
given by the Smoluchowski equation j^, 14 : 



d t ^! = fi^* j where ^J2d i -(d l -F i -Kr l ) . (1) 



The system is taken to be in quiescent equilibrium 
(7 = 0) at t < when averages (. . .)(t ,=0 ) are the canoni- 
cal equilibrium ones. Then at t = + , the velocity profile 
is switched on instantaneously, so that the steady state 
distribution function \& s , which satisfies f2W\& s = 0, will 
be approached at long times, t — > 00. If ^ s was known 
the steady state average a — (a xy )^' of the (thermody- 
namic) shear stress [|l5| could be found. From this the 
shear viscosity would follow as 77(7) = + cr/7 (where 
the solvent viscosity is denoted 7700). But the rapid in- 
crease of 77(0) close to the glass transition suggests that 
\l/ s is sensitively dependent on shear rate, which makes a 
direct calculation of it difficult. 

A more robust way to approximate steady-state quan- 
tities comes from the insight that the growth of 77 at the 
glass transition arises by slowing down of structural re- 
laxations H , whose characteristic time t in the quiescent 
state defines a second, "dressed" Peclet (or Weissenberg 
I]) number, Pe = jr. This characterizes the influence 
of shear on structural relaxation and increases without 
bound at the glass transition, even while Peo <C 1. We ar- 
gue that the competition of structural rearrangement and 
shearing that arises when Pe> 1 ^Peo dominates the 
non-linear rheology of colloids near the glass transition. 
Therefore steady state quantities shall be determined by 
considering the structural relaxation under shearing and 
"integrating through the transient dynamics". Because 
Peo <C 1, we expect ordering or layering transitions to 
be absent |7j ; and as hydrodynamic interactions are pre- 
sumed to play a subordinate role during the structural 
relaxation || we neglect these too, focusing solely on the 
Brownian contribution to the transverse (shear) stress. 

With t the time passed since the start of shearing, the 
correlation function of fluctuations in variables / and g 
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separated by time t is given by: 

C fg (t,to) = (e n » )t °r (e n *%)> (0) , (2) 

where the backwards Smoluchowski operator il^ ' arises 
from partial integrations: — ^,Ad% + F^ + rj/« T ) • di, 
and t,t > 0. Steady state expectation values, (g)^ = 
Ci g (t,t — > oo), and variances, (f*gp^ — Cf g (t 
0. £ — > oo) then follow. 

An important property of the sheared system is trans- 
lational invariance |jlq| . It leads to spatially indepen- 
dent averages, or, in Fourier space at wavevector q, to 
(/q)^ — (fo) ^q.o- I n the two-time correlation func- 
tions of Eq. (||), it leads to a coupling of fluctuations 
of wavevector q with later fluctuations of the advected 
wavevector q(f ) = q + q« t, suggesting the definition: 
^/kSq^o) = NC7s;q(Mo)<5 q (t),k- Figure H sketches 
the advection of a fluctuation with initial wavelength 
X x to one with X x and X y (t) — X x /(jt) at later time 
t. Brownian particle motions (assisted by the interaction 
forces) "smear out" the fluctuation with time and cause 
the decay of the corresponding correlator. Because the 
wavenumber q(t) increases upon shearing, smaller and 
smaller motions can cause the fluctuation to decay 

Equation (J2j) is an exact consequence of Eq. (Q) (for 
the given shear flow); but to proceed further requires 
additional approximations. With our assumption that 
applied shear interacts with slow structural rearrange- 
ments, we build on the description achieved by the MCT, 
and analyse the approach into the steady state by moni- 
toring the fluctuations of density (g q = J2iLi e jqri ) and 
of the "pair density" (the square of the density in real 
space) , aiming to establish nonlinear closed equations for 
these. This entails elimination of forces Fj in favor of the 
quiescent-state structure factor S q (taken to be known) 
- a near-equilibrium assumption (see below) that is for- 
mally uncontrolled but motivated, at least in part, by the 
smallness of Peo- 

Steady state correlators are now approximated by pro- 
jection onto the density modes, giving 

(3) 

with t now the time since switch-on, Si — dS^/dk, and 
V^ s an easily found static overlap fri. The transient 
density fluctuations are given by $ q (i) = C ee -q(t,0)/ S q , 
and are normalized by S q . They enter Eq. (|3|) via a 
factorization-approximation of the density pair fluctua- 
tion functions. With the choice / = g = gq/v^/V, Eq. 
(|J) gives the steady state structure factor under shear, 
whereas / = 1 and g — o~ xy /V give the transverse stress, 
for which Vl 3 = (Nk x k y /Vk)S' k . 

The problem of calculating the steady state averages is 
thus "reduced" to first finding the transient density fluc- 
tuations $ q (i), given by the structural rearrangements 



after switching on the flow, and integrating these in Eq. 
(0). From Zwanzig-Mori (type) manipulations [jis] , one 
finds the exact equation of motion: 

$ q (t) + r q (t) |$ q (t) + Jdt' m q (t, t') $ q (i')j = . (4) 

where < f > q (t) = <9 t $ q (£), and the "initial decay rate" r q (t) 
exhibits the familiar Taylor dispersion Jf4|, . It is not 
known how to evaluate the microscopic expression for the 
memory function m q (t,i') exactly. In the MCT spirit 
of our approach, this is approximated by projecting the 
fluctuating forces onto density pairs and factorizing the 
resulting pair-density correlation functions as 

m ^ 1 ') ra ^k(M') $k(*-t')* q -k(i-0 • (5) 

k 

The vertex , whose lengthy formula will be published 
elsewhere, is evaluated in the limit Peo <C 1 (as argued 
above) but for large times so that jt and jt' are finite. 
As 7 — > 0, it reduces to the standard MCT vertex |5) and 
like the latter is uniquely determined by the equilibrium 
structure factor, S q . For long times, it vanishes as " oc 

qyq-^-H'^a-t')- 1 forfe^O. 

Equations (g) to (||) complete our derivation of closed, 
self-consistent equations for the steady state properties 
of dense sheared suspensions. They contain the bifur- 
cation singularities which lie at the core of MCT. For 
7 = 0, upon smooth changes of the input equilibrium 
state parameters, a fluid with $ 9 (i — > oo) — > turns into 
an amorphous solid, $ g (t — > oo) — > f q > 0. The f q are 
called glass form factors and describe the arrested struc- 
ture. While transport coefficients of the fluid, like the 
viscosity, are connected to the longest relaxation time of 
& q (t), elastic constants of the solid, like the transverse 
elastic modulus Goo, are given by the f q g. 

In the limit of small shear rates, a stability analysis 
of the amorphous solid can be performed and leads to a 
generalization of the factorization theorem of MCT |20] . 
Close to the bifurcation, the dynamics are governed by 
®q(t) = fq + h q G{t), where the f q describe the glassy 
structure at the instability and the critical amplitude h q 
is connected to the cage-breaking particle motion; both 
retain their definition from the unsheared situation. Here 
Q{t) contains the essential non-linearities of the bifur- 
cation dynamics which arise from the physical feedback 
mechanism (the cage effect). It depends on a few mate- 
rial parameters only and, for \G(t) \ -C 1, follows from: 

e - c « ( 7 tf + xg 2 (t) = ^- t j dt' g(t ~ f ) g(t') . (6) 

Here, e measures the distance to the transition and A is 
known for some systems 0]. Monte Carlo estimations of 
the microscopic expression for the (new) parameter c^ 
give c^ 1 ' ~ 3 for hard spheres. Corrections of higher 
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order in the small quantities (e, 7, Q) are neglected; see 
PH f° r background on Eq. @ for 7 = 0. 

Equation (^J) is our central result. As expected, the 
sign of 7 does not enter, although it affects the Taylor dis- 
persion. Because {^t) 2 dominates for long times, always 
G{t->oo) -» -i/rW, with = v/(A- l/2)/c("V)/|7|. 
Hence, under flow, density fluctuations always decay, as 
this decrease of Q(t) for long times initiates the final re- 
laxation (where the corrections to Eq. (^) become im- 
portant) of $ q (£) to zero. Arbitrarily small shear rates 7 
melt the glass and so "interrupt" aging, as has also been 
found for spin-glasses where shearing was mimicked by 
breaking detailed balance O]. This vindicates our deci- 
sion to ignore aging and to proceed via Eq. (j^) in order 
to obtain steady state properties. 

While the non-Newtonian fluid behavior (e < 0) in- 
cludes two slow time scales, the familiar r and the shear 
induced the rheology of glass (e > 0) is determined 
by t^- 1 only ]22||. For e — > 0+ and 7 = 0, a fraction 
f q = fq+h q y/e/{\ — A) of the density fluctuations would 
stay arrested, while with shear these decay at a rate set 
by 7: $ q (t) -> /, $+(t/rW), where $+(x^0) -k-i, 
as follows from Eq. (||). 

For £ > and 7^0, the time for the final decay 
can become arbitrarily slow compared to the time char- 
acterizing the decay onto f q . Inserting <&+(t/M") into 
Eq. (||), the long time contributions separate out, and 
depend on time only via jt. Hence the glass has nontriv- 
ial shear-rate-independent limits for steady state values: 
(/*.9> W - (/*S> (+) ^ (/*.9> (0) for 7 -» 0. The (/*«,)(+) 
are given by integrals over the $q(t/r^- ) ) and quantify 
those fluctuations that require the presence of shearing 
to avoid their arrest. For the case of the shear stress, a 
finite (dynamical) yield stress, er + = ]hn^o{a xy )^ for 
e > 0, is thereby found. Since the glass transition is often 
identified by a divergence of viscosity (in terms of which 
we have shear thinning: 77(7 — > 0) oc -y _m with m = 1) 
our prediction of a finite yield stress throughout the glass 
phase is far from trivial. It excludes e.g. power-law fluid 
behavior (a ~ A f 1 ~ m with < m < 1; see [[To||). 

These results follow from the general stability analy- 
sis of Eq. (0) and are predicted to be universal, i.e. to 
hold for the Brownian contribution to the shear stress 
and viscosity close to the glass transition in any colloidal 
dispersion. But because (in contrast to aging approaches 
p3[ ) we approximate nonlinear couplings under shear us- 
ing equilibrium averages, we require the system to remain 
"close to equilibrium" in some sense. The existence of a 
finite yield stress cr + means that this is not guaranteed 
even as 7 — > 0. 

Insight into an important mechanism of shcar- 
fluidization can be gained by considering transient den- 
sity fluctuations with wavevector q perpendicular to the 
flow plane, q = qe z . Here simplifies to the standard 
MCT vertex with advected (time-dependent) wavevec- 
tors p3. While for 7 = it exhibits the nonlinear cou- 



pling of density correlators with wavelength equal to the 
average particle distance, for 7 7^ 0, the (only) effect of 
shearing on this "neutral" direction consists in a shift of 
the advected wavevectors to higher values, where the ef- 
fective potential decreases. This decreases the memory 
function and thus speeds up structural rearrangements. 
In this way the theory captures the faster decay of fluc- 
tuations caused by shear advection (cf Fig. |l]) . 

The presence of shear advection in the neutral z- 
direction suggests an approximation that considers only 
the resulting competition of caging and advection- 
induced decay. In this "isotropically sheared hard sphere 
model" (ISHSM) |25|] , we neglect kinematic flow of par- 
ticles so that all directions are treated as "neutral" . The 
quiescent S q depends only on the packing fraction 0, and 
the model's glass transition lies at <f) c = 0.51591 pi) . 
Figure || shows the stress versus strain rate curves for <fi 
close to the transition. In the fluid, <j> < <j) c , a Newto- 
nian regime {a = 777) is found for 7 small enough that 
Pe < 1. For Pe > 1 there is a broad crossover to the crit- 
ical yield stress value, = <j + (e — > 0+), from which a 
starts to rise due to (non-universal) short-time effects for 
Peo around 10 -2 (where t]^ will also contribute). In the 
glass, <f> > <f> c , a yield stress plateau for 7 — > is obtained, 
and rises strongly with increasing packing fraction. Wc 
speculate that the lack of a clear yield stress plateau at 
Peo > 1CP 3 (and likely hydrodynamic effects) explains 
shear-thinning exponents m < 1 seen in experiments |Q . 

In summary, we have presented a microscopic theory 
of the nonlinear rheology of colloidal fluids and glasses 
under steady shear. It predicts a universal transition 
between shear-thinning fluid flow, with diverging viscos- 
ity upon increasing the interactions, and solid yielding, 
with a yield stress that is finite at (and beyond) the glass 
point. Besides its interest for dispersion flow, our work 
suggests a further role of colloidal systems in elucidat- 
ing glasses via the study of shear-melted states for small 
shear rates. This is of fundamental interest because, e.g., 
mean-field driven spin-glass theories predict nothing like 
a yield stress p|. While comparison with measurements 
in colloids f| [7 , ||] and simulations of sheared atomic 
glasses Jl^, []L3) is promising, our approach represents only 
the first step toward rational prediction of the rheology 
of a glass; though physically motivated and in part in- 
spired by the successful MCT description of the cage- 
effect, several of our approximations remain incompletely 
justified. Nor is it clear whether Eqs. (§,||) can exhibit 
"jamming" transitions p6| at finite shear rates, or how 
anisotropic J2^J the fluctuations can become. Extension 
to time dependent shearing would be especially interest- 
ing because recent shear echo measurements [^8| reveal 
intriguing glass-melting scenarios. 
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FIG. 1: Advection by steady shear of a fluctuation in x- 
direction with wavelength X x at t = 0. At later time t, its 
wavelength X y in y-direction obeys: X x /X y = Ax/ Ay = jt. 
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FIG. 2: Steady shear stress, a = (o' xy }^\ in units of feT/d 3 
versus Peo = jd 2 /Do, for a model of hard spheres |25| at 
various distances from its glass transition, {(j> — (j> c ) as labeled. 
For the fluid cases, <f> < cj> c , dashed lines indicate Newtonian 
fluid behavior, a — 777, while vertical bars mark Pe= 7T = 1, 
with the structural relaxation time taken from ( I > ( j = 7/d(i = 
r) = 0.1. For the critical density, <f> c , the critical yield stress, 
crj" = 6.0, is shown by a horizontal bar, and the dotted line 
a = ct+(1 + I.O7 ' 17 ) matches for 7^0. 



